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Abstract. Numerical analysis of time-integration algorithms has been applying ad- 
vanced algebraic techniques for more than fourty years. An explicit description of 
the group of characters in the Butcher-Connes-Kreimer Hopf algebra first appeared 
in Butcher's work on composition of integration methods in 1972. In more recent years, 
the analysis of structure preserving algorithms, geometric integration techniques and in- 
tegration algorithms on manifolds have motivated the incorporation of other algebraic 
structures in numerical analysis. In this paper we will survey algebraic structures that 
have found applications within these areas. This includes pre-Lie structures for the ge- 
ometry of flat and torsion free connections appearing in the analysis of numerical flows 
on vector spaces. The much more recent post-Lie and D-algebras appear in the analysis 
of flows on manifolds with flat connections with constant torsion. Dynkin and Eule- 
rian idempotents appear in the analysis of non-autonomous flows and in backward error 
analysis. Non-commutative Bell polynomials and a non-commutative Faa di Bruno Hopf 
algebra are other examples of structures appearing naturally in the numerical analysis of 
integration on manifolds. 
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1 Introduction 

We present some aspects of the field Geometric Numerical Integration (GNI), 
whose aim it is to construct and study discrete numerical approximation methods 
for continous dynamical systems, exactly preserving some underlying geometric 
properties of the system. Many systems have conserved quantities, e.g. the energy 
in a conservative mechanical system or the symplectic structures of Hamiltonian 
systems, and numerical methods that take this into account are often superior to 
those constructed with the more classical goal of achieving high order. 

An important tool in the study of numerical methods is Butcher series (B- 
series), invented by John Butcher in the 1960s. These are formal series expansions 
indexed by rooted trees, and have been used extensively for order theory and 
the study of structure preservation. We will put particular emphasis on B-series 
and their generalization to methods for dynamical systems evolving on manifolds, 
called Lie-Butcher series (LB-series). 

It has become apparent that algebra and combinatorics can bring a lot of insight 
into this study. Many of the methods and concepts are inherently algebraic or 
combinatoric, and the tools developed in these fields can often be used to great 
effect. Several examples of this will be discussed throughout. 



2 Numerical integration on vector spaces 

In numerical analysis one of the main objects of study is flows of vector fields, 
given by initial value problems of the typcQ 

y'{t)=F{y{t)), y(0) = 2/ . (1) 

The function y can be real-valued or vector-valued (giving rise to a system of 
coupled differential equations). The flow of the differential equation is the map 
q t<F : R n -> R™ defined by y(t) = * t , F (y ) Q Note that F(y) = al/dt\ t= ^ t ^ F {y ). 
In many practical settings the vector field F is Hamiltonian, and their flows have 
several interesting geometric properties. We seek to construct good approxima- 
tions to flows, where 'good' can mean several different things, depending on the 
context. Sometimes what we want are approximations of high order, other times 
we need them to preserve some qualitative or geometric structure of the underlying 
dynamical system. Preserving geometric structure is particularly important when 
studying systems over long time intervals. An early illustration of this fact was 
made by Wisdom and Holman in |87j . where they computed the solar system's 
evolution over a billion-year time period using a symplectic method, making an 
energy-error of only 2 x 10 -11 . 

As there are several excellent introductions to geometric numerical integration 
on vector spaces we will not go into a detailed study here, but merely describe some 

1 Non-autonomous differential equations can also be written on this form by adding a compo- 
nent to the y vector 

2 Here we assume Lipschitz continuity of F for the flow to exist and be unique. 
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of the main ideas. The book [42] is the standard reference; other introductions 
can be found in [SH IS1 E3 ISD ISH] • 

The focus of this paper will be on some of the algebraic and combinatorial tools 
of geometric numerical integration, with particular emphasis on the tools we will 
utilize when studying flows on more general manifolds in Section [3| Lately, there 
has been quite a lot of interest in these algebraic aspects of geometric integration, 
and this has resulted in both an increased understanding of the field, and also of 
its relations to other areas of mathematics. 



2.1 Numerical methods and structure-preservation 

Consider an initial value problem of the form (JlJ : 

y'(t)=F(y(t)), y(Q)=y (1) 

representing the flow of the (sufficiently smooth) vector field F. A numerical 
method for (JlJ generates approximations y%, yi, 2/3, • ■ • to the solution y(t) at var- 
ious values of t. One of the simplest methods is the (explicit) Euler method. It 
computes approximations y n to the values y(nh), where n £ N and h is the step 
size, using the rule: 

y n +i=yn + hF(y n ). (2) 

This generates a numerical flow $^ approximating the exact flow \t of F . The 
accuracy of the method can be measured by its order: we say that a one-step 
method y n+1 = $>h(Vn) has order n if \&h(y) ~ ^h(y)\ = 0(h n+1 ) as h -> 0. 
Another way to put this is in terms of the curve traced out by the numerical flow: 
by comparing its Taylor series to the Taylor series for the curve of the exact flow 
term by term, we can read off the order of the method. The Taylor series for y 
has the form 

y{h) =y + hF(y Q ) + hi 2 F< '(y )F(y ) + 0(^ 3 ), (3) 
and we note that the Euler method is of order 1. 

Runge Kutta methods. The Euler method is an example of a Runge-Kutta 
method, a class of methods that are extremely prevalent in applications [321 [TT] . 
A Runge-Kutta method is a one-step method computing an approximation y\ to 
y(h) with j/o as input, as follows: 

Definition 2.1. An s-stage Runge-Kutta method for solving the initial value 
problem (JlJ is a one-step method given by 

s 

Y l =y + h^a lJ F{Y :j ), i = l,...s 

3 : x (4) 

2/i =y + h^hF{Y l ), 
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where bi,a,ij G R, h is the step size and s € N denotes the number of stages, i.e. 
the number of evaluations of F. 

A Runge-Kutta method can be presented as a Butcher tableau, which charac- 
terizes the method completely: 



Cl 


0-11 ■ 


■ a ls 


Cs 


a s i . 






h . 


■ b s 



The coefficients Cj = X)j=i a ij> apP ear explicitly only in integration methods for 
non-autonomous equations. The method is called explicit if the matrix {ctij} is 
strictly lower triangular, and diagonally implicit if it is lower triangular. Otherwise 
the method is called implicit. 

Example 2.2. We note that the Euler method HQ is the Runge-Kutta method 
with Butcher tableau: 
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Another well-known example is the explicit midpoint method: 

Vn+i =Vn + hF (y n + ^hF(y n ) \ , (5) 

given by: 












1/2 


1/2 
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An implicit method with good structure preserving properties is the implicit mid- 
point rule 

Vn + Vn+i . h ( y n + y n+ i 
2 =Vn+ 2 F { 2 

represented by the tableau 

1 l 

2 2 
1 

Given any number m, there exist Runge-Kutta methods of order m [TT]. Ver- 
ifying this involves expanding the methods into series involving the derivatives of 
F, and already at low orders the expressions get quite complicated. However, in 
Section |2~2"1 we shall see that the Runge-Kutta methods are special cases of Butcher 
series methods, and that one can find nice descriptions of the order theory and also 
structure preservation properties for numerical methods within this framework. 

Differential equations and geometric structures. When presented with a 
system modeled by a differential equation one will often first try to determine its 
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qualitative properties: are there any invariants? What kind of geometric structure 
does the system have? Structures of interest can be energy and volume preserva- 
tion, symplectic structure, first integrals, restriction to a particular manifold (as 
studied in Section [3]), etc. Then, when choosing (or designing) a numerical method 
for approximating the solution of the differential equation, it might make sense for 
the method to share these qualitative features. In that way one has control over 
what kind of errors the method introduces, obtaining a method tailor-made to the 
problem at hand. 

A rich source of problems with geometric structures are the Hamiltonian sys- 
tems. Let H : M. 2n — > K be a smooth function. A Hamiltonian vector field is 
a vector field on E 2 " of the form Xjj = f2 _1 V.ff , where fi is an antisymmetric, 
invertible 2n x 2n matrix^ The flow of Xh is given by 



The function H represents the total energy of the system. Two important prop- 
erties of the flow of a Hamiltonian vector field Xh is that it is constant along 
the Hamiltonian function H (conservation of energy) and that it preserves a sym- 
plectic form ui on R 2n . Using numerical integrators constructed to preserve these 
properties has been shown to lead to dramatic improvements in accuracy. For 
examples of this phenomenon see e.g. |42l |41j [53] , and references therein. 



Starting with the work of John Butcher in the 1960s and 70s [HI HH] the study 
of methods for solving ordinary differential equations has been closely connected 
to the combinatorics of rooted trees. Many numerical methods y n +i = ^h{yn) 
(including all Runge-Kutta methods) can be expressed as certain formal series, 
called Butcher series by Hairer and Wanner in }44] , By a clever representation 
of the terms the series can be indexed over the set of rooted trees. 
Consider the differential equation 



^-z = nV z H(z). 



(7) 



2.2 Trees and Butcher series 



y'(x) = F(y(x)). 



(8) 



Denote the components of F : R" — > R" by f and write 



Qkfi 



(9) 



dxj 1 dxj 2 ■ ■ ■ dxj k 



3 Hamiltonian vector fields can be defined on any symplectic manifold [1]. 
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Summing over repeated indices, the first few derivatives of y can be written as: 
f 

fir 



dx 



dx 2 

fififtf + i]fij k f + zf) k fif k f + r jk if j f k f 



d V _fi j-j rk , ft rj fk 



d 4 y 



(10) 



dx 4 

These expressions soon get very complicated, but the structure can be made much 
more transparent by observing that the derivatives of F can be associated in a 
bijective way with rooted trees, an observation already made by Cayley in 1857 
|18j . Before giving the exact correspondence between differential equations, rooted 
trees and Butcher series, we will take a closer look at trees. 

Rooted trees. A tree is a connected graph with no cycles 

A rooted tree is a tree with one vertex designated as the root. In the pictorial 
representation of trees, the root will always be drawn as the bottom vertex, and 
the trees will be ordered from the root to the top. More precisely a tree r is a 
graph consisting of a set of vertices V(t) and edges E(t) C V(t) x V(t) so that 
there is exactly one path connecting any two vertices. A path between Vi and Vj 
is a set of edges {v Sl , i?t, } so that I = 1, 2, . . . , r, s\ = i, ti = s; + i and t r = j. This 
gives a partial ordering of the tree in terms of paths from the root to the vertices 
of the tree. A vertex V{ is smaller than another distinct vertex Vj, e.g. Vi -< Vj, if 
the unique path from from the root to Vj goes via Vi. A vertex Vi is called a leaf 
if there is no vertex Vj with Vi -< Vj. A child of a vertex vi is a vertex Vj with 
Vi -< Vj so that there is no vertex Vk with Vi -< Vk -< Vj. The order \t\ of a tree 
r is the number of vertices of the tree. We define a symmetry group on a tree r 
as all automorphisms on the vertices. The order of this group, <t(t), is called the 
symmetry of the tree r. 

A forest of rooted trees is a graph whose connected components are rooted 
trees, e.g. w = T\...T n . We include the empty tree I, i.e. the graph with no 
vertices, in the set F of forests. F can be put in bijection to the set of trees via 
the operator B + : F — > T, defined on a forest uj = t\ . . . T n by connecting the trees 
to a new root by addition of edges. For example, 



B 



This operator can be used to generate all trees recursively from the tree • by the 
following procedure: 

(i) The graph • belongs to T 
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(ii) If ti, . . . , T n € T then t = £> + (ti . . . r n ) is in T. 
The tree factorial r! is given recursively by: 

(i) .! = 1 

(ii) B+( n . . . T B )! = \B+(n . . . r„)|n! . . . r n !. 

An important operation on trees is the Butcher product, defined in terms of graft- 
ing. 

Definition 2.3. The Butcher product tooj of a tree r = B + (ti . . . r„) and a forest 
a; = u>i ■ ■ ■ u) m is given by grafting u) onto the root of t: 

TOW = B + (Tl . . ,T n (jJl . . .UJ m ) (If) 

Butcher series. The calculations of the derivatives of y'(t) = F(y(t)) performed 
at the beginning of the section can be written in terms of the elementary differen- 
tials of F. 



Definition 2.4. Let F 

T of F is 



•F(.)(t) = F(y) 
J-(r)(t) = FM (y ) (7 - (n )( 2/ ) 



" be a vector field. The elementary differential 

(12) 



where i^™) is the m-th derivative of the vector field F and r = £? + (ti, . . . , t to ) is 
a rooted tree. 



We will discuss another way to write elementary differentials in Section [275] With 
the notation from Equation the first few elementary differentials are shown 
in Table (1.1). The vector field F corresponds to the leaves of the tree, the first 
derivative F' corresponds to a vertex with an edge with one child, the second 
derivative F" corresponds to a vertex with two children, etc. 



T 




• 


r 


I 


f)f 3 


Y 


r jk pf k 


I 


rjir 




f l jk i.Pf k f l 




p jk f j ftf l 



Table 1. Elementary differentials associated to a vector field F with components f 1 . 
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Butcher series are (formal) Taylor expansions of elementary differentials in- 
dexed over trees: 

Definition 2.5. A Butcher series (B-series) is a (formal) series expansion in a 
parameter h: 

,a(r) 



B h>F {a) = a{l)F(l) + V h^^^r) 



(13) 



where T = TU{I}, F is a vector field, a is a function a : T — > R, <t(t) is the 
symmetry of r, /i is a real number (representing the step size), and J 7 is the 
elementary differential of F, extended to the empty tree I by J-"(I)(y) = y. 

We shall see that these series can be used to represent numerical methods 
y n +i — &h(yn) approximating the flow of a vector field F, in the sense that the 
Taylor series for 3> h can be expanded into a B-series: <&h — Bh.F{a)^\ 

By computing the Taylor expansion of the solution to the initial value problem 
([!]) one obtains the following result: 

Proposition 2.6 ( 42J). The Taylor series for the solution of the differential equa- 
tion |I]) can be written as a B-series: 

B h Ai) = Y, hlTll r\^> ( 14 ) 

where 7(7-) = 1/t!. That is, y(t + h) = Z^i? (7) (?/(£)). 

Runge-Kutta methods can also be written as B-series expansions, with coeffi- 
cients given by the elementary weights of the method [S] . 

Definition 2.7 (Elementary weights). Let hi and a,ij be coefficients of a RK- 
method as in Definition |2.1[ where i £ N. The elementary weight function $ is 
defined on trees as follows: 

$*(•) = Ci 



3 

3=1 



MB + (n, . . . ,r k )) =^ay$ J -(r 1 )^(T 2 ) . . .^(rfe) 



(15) 



*{B + (n,. ..,r k ))=J2 ^^(n)^(r 2 ) • • • %(r fe ) 



4 A numerical method for solving a differential equation is called a B-series method if it can 
be written as a B-series. 
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Here i = 1, . . . , s. 



For example, 

= J2 *(V) = ^ & jC f , $(Y) = ^ 6 jaj -*cg 

Theorem 2.8 ( 9 ). XTie B-series for a RK-method given by the elementary weights 
$(r) is 

^f(1>) = $> M$ ^(t) (16) 



Order theory for B-series methods. Once we have the B-series of the exact 
solution and the B-series of a numerical method, it is straightforward to compare 
the coefficients and read off the order of the method. For Runge-Kutta methods, 
we obtain the following result: 

Proposition 2.9 (,9s)- A Runge-Kutta method given by a B-series with coeffi- 
cients $(t) has order n if and only if 

$( T )= 7 (r) ; for all t E T such that |r| < n. (17) 



B-series methods and structure preservation. The class of B-series meth- 
ods includes all Taylor series methods and Runge-Kutta methods. It does not, 
however, include all numerical methods, an example being the class of splitting 
methods. 

It is important to point out that focusing only on B-series methods has its 
drawbacks. Besides the fact that the class does not contain all methods, it is also 
known that there are certain geometric structures that cannot be preserved by 
B-serics methods. For example, no B-series method can preserve the volume for 
all systems [SU]. However, we will be content with this loss of generality and focus 
exclusively on methods based on B-series in this section, and on their generalization 
- Lie-Butcher series - in the next. 

A case which is particularly well-studied is Hamiltonian vector fields. The 
following two theorems serve as prime examples: 

Theorem 2.10 ([30]). Let G = B, hF (a) be a vector field with a (I) = 0, a(») ^ 0. 

Then G is Hamiltonian for all Hamiltonian vector fields F(y) = il _1 ViJ(y) if and 
only if 

a(riOr 2 ) + afoon) = (18) 



for all Ti,T2 G T. Here o denotes the Butcher product of Definition 2.3 



Theorem 2.11 (|16j). Consider a numerical method given by a B-series Bh,F{o)- 
The method is symplectic if and only if 

a(n o r 2 ) + a(r 2 o n) = a(n)a(r2) (19) 
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for all Ti,T2 G T, where a(T) = 0. 

The paper |2()j gives an overview of what is known about structure preservation 
for B-series, including characterizations of the various subsets of trees correspond- 
ing to energy-preserving, Hamiltonian and symplectic B-series. 

2.3 Hopf algebras and the composition of Butcher series 

Consider two numerical methods given by and $ 2 . Using the method $ 1 to 
advance a point yo to a point yi, and then applying the method $ 2 using y\ as 
initial point, results in a point yi\ 

yi = $ 1 (yo), 2/2 = $ 2 (2/i)- 

This is the idea behind composition of numerical methods. In the case where both 
methods are given by B-series, ^(yi) = Bj lF (a)(y Q ), $ 2 (yi) = Bl F ((3)(y ), the 
composition method ^oQ 1 is again a B-series: $ 2 o$ 1 (y ) = Bh,F(l)(yo)- This is 
the Hairer- Wanner theorem from [33] . The coefficient function 7 of this B-series 
was first studied by John Butcher in [10] . where he found that composition of 
B-series is a group operation (giving rise to the Butcher group) on the coefficient 
functions, and gave expressions for the product, identity and inverse in this group. 

In [52l ES] Connes and Kreimer introduced a Hopf algebra of rooted trees 
connected to the renormalization procedure in quantum field theory. Later [7] 
it was pointed out that a variant of this Hopf algebra is closely related to the 
Butcher group. More precisely, the Butcher group is the group of characters in a 
Hopf algebra Hbck defined by Connes and Kreimer. 

We will describe the Butcher group indirectly by describing the Hopf algebra 
Hbck- But first we will present some basic definitions from the theory of Hopf 
algebras. For a comprehensive introduction, see [5U1 El- Other excellent references 
include [T71 [BO]. 

Hopf algebras. Let k be a field of characteristic zero. An algebra A over k is 
a k-vector space equipped with a multiplication map /i : A ® A — > A and a unit 
u '. k — y A so that 

• /j, o [id ® /i) = /-i o (/i <S> id) : A® A® A — > A (associativity) 

• /1 o (u ® id) = n o (id <Eiu) : k ® A = A — >• A (unitality) 

A coalgebra C over k is the dual notion. It consists of a comultiplication map 
A : C — > C <g> C and a counit e : C — > k so that 

• (A (g) id) o A = (id ® A) o A : C -> C ® C <£> C (coassociativity) 

• (e ® id) o A = (id ® e) o A : C ->• C <g> k = C (counitality) 

A Hopf algebra is at once an algebra and a coalgebra, and it comes equipped with 
an antipode S : H — > . These structures have to satisfy certain compatibility 
conditions, written as the following diagrams, where r denotes the flip operation 
r(h ll h 2 ) = (h 2: hi): 
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H 



H 



H®H — £ k®k 



H 




1®S 



The first two diagrams ensure that the coproduct and the counit are both algebra 
homomorphisms. The last diagram is best interpreted in terms of the characters 
in a Hopf algebra. Let A be a commutative fc-algebra, and let C(H 1 A) denote the 
set of linear maps from H to A. An element a £ C{H, A) is called a character if 
a(x ■ y) — a(x) ■ a(y) for all x,y £ H, where the product on the left-hand side is 
in H, and on the right-hand side in A. The set of characters in C(H,A) form a 
group under the convolution product: 



(f) * 1p = ^ o (0 ® ip) o A. 



(20) 



The unit is the composition of the unit and the counit in H, e.g. r/ := u o e. 
The bottom diagram above corresponds to the antipode being the inverse of the 
identity under this product, and we have a* -1 = a o S. 

Later we will also need the concept of infinitesimal characters (also called 
derivations), which are maps a in C(H,A) satisfying 



a(x ■ y) = i](x) ■ a(y) + a(x) ■ r)(y). 



(21) 



The Butcher Connes Kreimer Hopf algebra. The composition of B-series 
is governed by a certain Hopf algebra Hbck based on the set T of rooted trees, 
called the Butcher- Connes- Kreimer Hopf algebra. In Section [3] we will see that a 
generalization of this Hopf algebra governs the composition of Lie-Butcher series 
(Section 2|3li|). 

To describe the BCK Hopf algebra we need to define its structure as a vector 
space, an algebra, a coalgebra, and define the antipode. As a R-vector space Hbck 
is generated by the set T of rooted trees, and graded by the order (i.e. number 
of vertices) of the trees. The algebra structure is that of the symmetric algebra 
^(^{T 1 }). The product is written as (commutative) concatenation of trees (i.e. 
disjoint union), giving rise to forests of trees. The unit is the empty tree I. 

iv=v: i i=i i=i 
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The coproduct of Hbck is the map Abck : Hbck — > Hbck ® Hbck determined 
recursively by: 

Abck ° B + M = ® I + (Id ® B+) o A BC kM, (22) 

where w is a foreslQ The counit is the map e : Hbck — > K given by e(I) = 1 and 
e(r) = if t =/= I. The coproduct can also be written in a non-recursive manner 
using cuttings of trees. 

Cutting trees. An admissible cut of a tree r is a set c C E(t) of edges of r 
such that c contains at most one edge from any path from the root to a leaf. The 
case c = is called the empty cut. The cut c = E(t) is allowed, and is called the 
full cut. Let uj denote the forest with vertices V(r) and edges E(t) \ c. We write 
R c (t) for the component of U) containing the root of r, and P c {t) for the forest 
consisting of the remaining components. 

Theorem 2.12 ('|26j). The coproduct in Hbck can be written as 

A BC k(t)= Y, P c (t)®R c (t) (23) 

c^Adm(r) 

Examples of the coproduct can be found in Table [2j The antipode can be defined 



Abck(t) 



I 
I 

V 



Y 



Y, 



+ »«(g)« + 2»(g)I + I(g)V 



■•® V + 2«® I + M( 



il+i® Y 

il + l® ¥ 



Table 2. Examples of the coproduct Abck in the Hopf algebra Hbck 



5 Recall that Abck is an algebra morphism and is therefore defined on forests as well as trees, 
since AbckOits) = ^bck( t i)^bck( t 2)- 
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recursively as: 



S(r) = -r- S(P c (t))R c (t) (24) 

cGAdm(r)\0 



The Hairer-Wanner theorem gives the exact correspondence between Hbck 
and composition of B-series: 

Theorem 2.13 (|44|). Let Bl F (a) and B F (f3) be two B-series, with coefficients 
a,P : T — > M. The composition B\ F {j3)oB\ F {a) is again a B-series, and we have 

Bl F (f3)oBl F (a) = B h , F (a*P), (25) 
where * denotes convolution in the Hopf algebra Hbck- 



2.4 Substitution and backward error analysis for Butcher 

series 

Consider a numerical method <&h used to solve a differential equation of the form 

y' = F(y). (26) 

The basic idea of backward error analysis of the method &h is to interpret it as 
giving the exact solution of a modified equation: 

y' = F h {y). (27) 

If we can find such an equation, we can use it to study the properties of the 
numerical method. In other words, the numerical method <&h will be represented 
by a modified vector field F, which then can be used to study the method. The idea 
is based on work by Wilkinson in the context of algorithms for solving equations 
given by matrices [55] , and has been explored in several papers [5SJ 001 E3 EH [H] . 
Recurrence formulas for the modified equation was first obtained in [ICJ [T5] . 

A related notion is the modifying integrators of |25j . The idea is to look for a 
vector field Fh so that the numerical method applied to the flow equation of 



20 



Fh (Equation 27 ) is the exact solution of Equation 

It turns out that the case where $/j is a B-series method is particularly nice 
[2"4"1 [231 [T2"] . The vector fields Fh can then be written as B-series whose coefficients 
are derived from the coefficients of $h, an d these coefficients can be expressed by 



the substitution law for B-series methods (Corollaries 2.15 and 2.10). 



The substitution law. Let Bh p(ce) and Bh cifi) De two B-series, where a(T) = 
0. Then Bh pip) is a vector field, and we can consider the B-series obtained by 
using this as the vector field G in the B-series Bh,G(P)- This is called substitution 
of B-series. The result is given in terms of a bialgebra Hcefm by the following 
theorem: 



Theorem 2.14 ([H]). Let F be a vector field, a, linear maps a, (3 : T — >• R 
where j3 is an infinitesimal character of Hbck, and a(I) = 0. Then the vector 
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field (l/h)Bh,F{&) inserted into the B-series Bh,-{/3) is again a B-series, given by 
Bh,(i/h)B h , F ( a) (P) = B hiF (a * (3), (28) 
where * denotes convolution of characters in the bialgebra Hcefm- 

The bialgebra Hcefm is the symmetric algebra over rooted trees S(T), with • 
as unit, equipped with a coproduct given by contracting subforests in trees: 



A(r) 



t/uj. 



(29) 



If t is a tree then the notation co C r means that w is a spanning subforest of r, i.e. 
that co is a collection of subtrees of r so that each vertex of r belongs to exactly 
one tree in co. Then t/uj denotes the tree obtained by contracting each subtree 
(with at least two vertices) of r contained in to onto a vertex. Some examples of 
the coproduct can be found in Table [3j 

There is a Hopf algebra related to Hcefm, obtained by considering the sym- 
metric algebra over the set of rooted trees T with at least one edge (e.g. • is not 
included), and then adding • back as the unit for the product. The coproduct is 



defined as in Equation ( 29 1 . This makes the associated bialgebra connected, and 



it is therefore a Hopf algebra [60] . 

For details on these constructions, consult [T2] . 



Acefm(t) 



I 

V 



Y 



I <g) • + m ® I 

I (g> • + M« (g) l 



I (g)»- 



Y + 2l»®I 



1. 



+ !>:'>!•• i n I 

Y+2iw(g)i+t^(g)V+2L«) 

^/ + I^(g)I + 2t^(g)V+ V« ® ! + L ® 2 



Table 3. Examples of the coproduct Acefm in the substitution bialgebra 



Backward error analysis and modifying integrators. Once Theorem |2.14| 
is established one can obtain expressions for backward error analysis and modifying 
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integrators. 

Corollary 2.15 (Backward error analysis). Let denote the B-series for the 

exact flow of the vector field G, and let Bf(oi) be a B-series giving a numerical 
flow for F . The modified vector field F given by B F {^f) = Bp{a) is a B-series 
Bf(P) with coefficients given by 

(3*j = a (30) 

Corollary 2.16 (Modifying integrators). Let Bail) denote the B-series for the 
exact flow of the vector field G, and let Bf(ch) be a B-series giving a numerical 
flow for F. The modified vector field F so that Bp(a) = Bf{i) is a B-series Bf{P) 
whose coefficients are given by 

(3 * a — 7 (31) 



2.5 Pre-Lie Butcher series 

The space of vector fields on R n has the structure of a pre-Lie algebra, and in this 
section we will see that B-series can be formulated purely in terms of this pre-Lie 
structure. This allows us to lift the concept of B-series to the free pre-Lie algebra, 
giving rise to pre-Lie B-series |31j . Viewing B-series as objects in the free pre-Lie 
algebra gives a clearer focus on the core algebraic structures at play, and it also 
enables the application of tools and results from other fields where pre-Lie algebras 



appear. Examples of this phenomenon can be found in [30] (see Remark 2.231, [12] 
and [2JJ. We give the basic constructions here because formulating Butcher series 
in terms of pre-Lie algebras will find an analogue in Section [3] where Lie-Butcher 
series will be constructed from the so-called post-Lie algebras. 



Pre-Lie algebras. The concept of pre-Lie algebras is a relaxation of associative 
algebras that still preserve their Lie admissible property. In other words, for an 
associative algebra (A, *) antisymmetrization of the product * gives a Lie bracket, 
making it a Lie algebra: [a, b] — a*b — b*a, and this property also holds for pre-Lie 
algebras. Note, however, that not all pre-Lie algebras are associative. They were 
first introduced and studied by Vinberg [S3] , Gerstenhaber [3S| , and Agrachev and 
Gamkrelidze [3J, under various names. A nice introduction to pre-Lie algebras can 
be found in [ST]. 

Definition 2.17. A (left) pre-Lie algebra^ (A, >) is a fc-vector space A equipped 
with an operation > : A <E> A — > A subject to the following relation: 

a&(x,y,z)-a > (y,x,z) = 0, (32) 

where a > (x, y, z) is the associator a>(x, y, z) = x > (y > z) — [x > y) > z. 

Example 2.18 (The pre-Lie algebra of vector fields). The space of vector fields 
X(M) on a differentiate manifold M equipped with a flat, torsion-free connection 



6 Also called a Vinberg, left- symmetric or chronological algebra 
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V can be given the structure of a pre-Lie algebra by defining > as F > G = V pG. 
In the case M = K" with the standard fiat and torsion-free connection we have 
that for F = Yh=i f A and G = 2"=i G fiv 

Ti I 71 \ 

^>g=E E F ^' G<) r' (33) 

In Section^ we will see that allowing for torsion leads to the concept of post-Lie 
algebras. See also \6&j. 



The free pre-Lie algebra. The free pre-Lie algebra has been studied in sev- 
eral papers, most notably by Chapoton and Livernet in [23] . Dan Segal in |77j . 
Agrachev and Gramkrelidze in [3], Dzhumadil'daev and Lofwall in [55]. These 
papers give different bases for the free pre-Lie algebra, and one can choose to work 
in the basis most beneficial for the problem at hand. A basis for the free pre-Lie 
algebra PL(V) over a vector space V was described by Chapoton and Livernet in 
terms of nonplanar rooted trees [21] : 

(.:v:t%>Y j 

decorated by elements of V. The pre-Lie product t\ rv T2 of two rooted trees is 
given by grafting: t\ rx is the sum of all the trees resulting from the addition 
of an edge from the root of T\ to one of the vertices of t-i'- 

T X rxT 2 := ^2 n o v t 2 (34) 

v£V(t 2 ) 

Here t\ o v t 2 denotes grafting at the vertex v of t%. 

< i 

• rv« = I, «rvI = V + I, IrvI = 2V + !! 

Theorem 2.19 f |23j). PL(V) is the free pre-Lie algebra on the vector space V: 
for any pre-Lie algebra P equipped with a morphism V — > P, there is a unique 
pre-Lie morphism PL(V) — > P making the following diagram commute: 

V PL(V) 



P 

We write PL for the free pre-Lie algebra on a space with only one element. 
The free pre-Lie algebra is related to the Hopf algebra Heck defined in Section 

I2~3l 



Theorem 2.20 ([23]). The universal enveloping algebra U(PL) of the free pre-Lie 
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algebra on the one-vertex tree, viewed as a Lie algebra, is isomorphic to the dual 
of the Butcher-Connes-Kreimer Hopf algebra Hbck- 

In fact, the dual of the Butcher-Connes-Kreimer Hopf algebra is isomorphic 
to the Grossman-Larson Hopf algebra defined |39j . The isomorphism was proven 
in [45]. 

Pre-Lie Butcher series. Now we can formulate the pre-Lie Butcher series 

Definition 2.21. A pre-Lie Butcher series is a formal series in R(PL): 

X(a) = h W a(t)t. (35) 



The classical B-series are recovered by applying the unique pre-Lie morphism 
associated to a vector field F: 

T : PL -> X(R n ) such that = F. (36) 



This is the elementary differential function of F as defined in |2.4| It is given 
recursively by J r (») = F and 

Ht)=F^{T(T l ),...,T{T n )), (37) 



if < = B + (t 1; ...,t„). 

B-series in any other pre-Lie algebra {A, >) can be defined in the same way: by 



applying the unique pre-Lie algebra morphism F : PL — > A to the series (35 1. 



Remark 2.22. Since J- : PL — > X(W l ) is a pre-Lie morphism, the trees associated 
to the derivatives of y'(t) = F(y(t)) can be generated by iterated grafting onto the 
one- vertex tree: 

d n y 

• rv (• rv (• rv . . . (• rv •) . . . )) corresponds to 

This way of looking at elementary differentials will reappear in a different setting 
in Section |3] 



Remark 2.23. The formulation of differential equations in terms of pre-Lie alge- 
bras has seen some use in numerical analysis. In |30) Ebrahimi-Fard and Manchon 
rephrased differential equations of the type X'(t) = A(t)X(t), where X, A are 
linear operators in a vector space, as combinatorial equations in pre-Lie algebras. 
In this context they obtained an analogue of the Magnus expansion [59], a series 
expansion of the solution to the equation in the magma generated by monomials of 
pre-Lie elements. In this setting it becomes apparent that one can use the pre-Lie 
relation to cancel out some of the terms in the expansion, leading to a hitherto 
unknown reduction of the number of terms in the Magnus expansion 
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3 Geometric numerical integration on manifolds 

Our objects of study are now dynamical systems evolving on manifolds: 

y'=F(y), y eM, F e X(M), (1) 

where M is a smooth manifold and X(M) denotes the vector fields on M. As 
in the previous chapter, the aim is to find good numerical approximations to the 
flow exp(tF) :— &t,F of |l]). The study of such systems comprises several different 
approaches: One simple way to attack the problem is to embed the manifold in 
R , for some AT, and use methods developed for M. N to solve the equation. But 
then the numerical flow of the method may drift off the manifold, and this can in 
some cases cause problems [341 |46l [13l [51] . 

A more satisfying and often better way is to use methods that are intrinsic to 
the manifold, and not rely on any embedding, which is the approach taken by Lie 
group integrators |48) . Consider for instance a system evolving on the manifold 
S 3 . By embedding S 3 in R one can use numerical methods that approximate the 
flow of the system using the basic motions of translations in R . Another approach 
is to use rotations to move around S 3 : y n +i = QnUn where Q n are orthogonal 
matrices, i.e. to use the action of the Lie group SO (4) on S 3 . This illustrates 
the intrinsic approach, where we are guaranteed not to drift off S 3 . Methods 
developed for manifolds include the Crouch-Grossman and RKMK-methods (and 
variants thereof) |g51 E3 E2 E2 E2 ■ 

In this chapter we will study a generalization of B-series called Lie-Butcher 
series. In analogy to the previous chapter we will look at the composition and 
substitution of Lie-Butcher series. 



3.1 Setting the stage: homogeneous manifolds and 
differential equations 

The flows we would like to approximate evolve on smooth manifolds, and so the 
tools of differential geometry play an important role. We will not review the gen- 
eral theory of smooth manifolds here, but assume a basic knowledge of differential 
geometry; for excellent introductions see e.g. [751 [75]. For a viewpoint oriented 
toward geometric numerical integration, see [48] . More precisely, we will be work- 
ing with smooth manifolds equipped with transitive actions by Lie groups, so called 
homogenous manifolds, where the Lie group provides a way to move around on 
the manifold^ Because the action is not in general free, the differential equation 
expressed on the Lie group is not in general unique. 

Definition 3.1. An action of a Lie group G on a smooth manifold M is a group 
homomorphism A : G — > Diff(M), g H> X g , where Diff(M) is the group of diffco- 
morphisms on M. We will mostly write such an action as a map A : G X M — > M. 



7 Note that other manifolds with local actions could also be considered, but to to avoid un- 
necessary complications we only consider homogeneous manifolds. 
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For convenience of notation we write g for the diffcomorphism X g , and also g ■ m 
for X g (m). The orbit through a point p 6 M is the set G-p = Xq(p). The action is 
called transitive if the manifold M is a single G-orbit. That is, if for all p,q £ M 
there is a g £ G so that p — g ■ q. A manifold equipped with a transitive action 
by a Lie group G is called a homogeneous manifold. A consequence of this is that 
M is diffeomorphic to the right cosets G/G x of G, where G x is the closed Lie 
subgroup of isotropies, G x — {g € G \ gx = x} (the point stabilizer): the smooth 
manifold structure of G/G x comes from the quotient map, and the diffcomorphism 
F : G/G x — > M is given by F(gG x ) — g ■ x. The group G x is called the subgroup 
of isotropies because if x' is another point in G, then G x and G x > are conjugate, 
and therefore isomorphic. 

Some interesting examples of homogeneous manifolds are the spheres S n — 
SO(n + \) / SO{n). Other important examples come from Lie groups G them- 
selves. The action A : G x G —> G is then the Lie group multiplication. A 
(somewhat degenerate) example is the homogeneous manifold (M n , (R™, +)). Here 
the action of R" on itself is given by translations. The theory developed for ho- 
mogeneous manifolds in this chapter will reduce to the theory developed in the 
previous chapter when applied to this particular case. 

Actions by Lie groups on manifolds can be associated to actions by Lie algebras. 
Let A : G x M — > M be an action of G on M. The associated Lie algebra action 
A* : g — » X{M) of g on M is the Lie algebra anti-homomorphism defined by: 



A*(u)(p) = j f 



A(exp(^),p). (2) 

t=o 



This satisfies A*([J7, V]) = — [A*(i7), A*(V)], where the bracket on the right is 
the Jacobi bracket of vector fields. We sometimes write v ■ y for the element 
A*(u)(y) € T y M. The Lie-Palais theorem |73j ensures us that as long as the Lie 
group G is simply connected, then every action by g comes from an action by G|^] 
If F G X{M) is a vector field, then an element v so that X*(v) = F is called an 
infinitesimal generator for F. 



Remark 3.2. In some cases it makes sense to use other maps <fi : g — >• G (satisfying 
</>(0) = e and <f>'(0) = V) besides the exponential map to construct maps g — > 
X(M) as in Equation An overview of various maps of this kind, and their 
usefulness, can be found in [52] . 



Differential equations in homogeneous manifolds. Consider the differential 
equation on a homogeneous manifold {M. G, A): 

y' = F(y), j/qGM, F:M^TM. (3) 

The solution is the flow ^t,F — exp(tF) of the vector field F. The vector field can 
be written in terms of its infinitesimal generator as F — A* (v) : M — > TM for an 
element v £ g, and the transitivity of the action also allows us to construct a map 

s If the Lie group is not simply connected, then we can only lift the g-action to the universal 
covering group of G. 
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/ : M — > g so that 

F(y) = K(f(y))(y) = f(y)-y (4) 

Note that as long as the action is not free, this / is not unique: if / : M — > g is 
such a map, then / + i : M — > Q, where i(p) is in the isotropy subalgebra g p of 
g, is another map of the same type. This choice of isotropy class can be helpful 
when constructing numerical integrators [54] . 

The differential equation ^ can be written as: 

y' = fiv) • V, where f : M —> Q, (5) 

and this is the type of differential equation we will consider in this chapter. Note 
that in the classical case of (R 71 , (R n , +)), the exponential is exp(u) = v, and 
Equation [5] reduces to the ordinary differential equation Q . We also note that 
the class contains the equations formulated in terms of frames. 



Remark 3.3 (Frames and differential equations). In the literature for numerical 
integration of differential equations on manifolds the equations are often simplified 
by using a frame on the manifold [75J [TTJ [T5] . A frame is a set of vector fields 
{Ei} that at each point on the manifold spans the tangent space at that point, so 
that any vector field F can be written as F = 'J2 i fiEi. The flow equation ^ for 
F can then be written as 

y' = 2_.fi(y)Ei(y), where /, : M —> M are smooth. (6) 

i 

If we write g C X(M) for the Lie subalgebra generated by the frame vector fields 
{E^, and let A* : g — >• Diff(M) be as in ([2]), we see that Equation is a special 
case of Equation ([5|, with / : M — > g defined by f(y) = £\ fi(y)E{. 



Remark 3.4. In j32], K. Eng0 discusses the general operation of 'moving' differ- 
ential equations between manifolds using equivariance of actions and relatedness 
of vector fields. In particular, every differential equation of the form ^ is shown 
to be equivalent to a differential equation on g. The following diagram from [32] 
summarizes this: 



Tgl^lTG 







G 



T(A.(p)) 



MP) 



TM 



M 



In other words, the differential equation on a homogeneous manifold (M, G) is 
moved to the Lie group G (the middle vertical arrow) and then to the Lie algebra 
(the first vertical arrow). As before, the exponential map exp : g — > G can in 
many cases be replaced by other maps. The construction of the vertical arrows 
can be found in This is the result exploited in the so-called RKMK methods 

15511551153. 
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3.2 Lie group integrators in applications 

Definition 3.5. Given a smooth manifold M and a Lie group G with Lie algebra 
g acting on M. Consider a differential equation for y(t) £ M written in terms of 
the infinitesimal action as 

m = f&v)-v, y(o) = j/o, (7) 

for a given function / : ExM — > g. Note that we now make the non-autonomousness 
explicit. A Lie group integrator is a numerical time-stepping procedure for Q 
based on intrinsic Lie group operations, such as exponentials, commutators and 
the group action on M. 

Applications of LGI generally involve the following steps: 

(1) Choose a Lie group and Lie group action that can be computed fast and 
which captures some essential features of the problem to be solved. This is 
similar to the task of finding a preconditioner in iterative solution of linear 
algebraic equations. 

(2) Identify the Lie algebra, commutator and exponential map of the Lie group 
action. 

(3) Write the differential equation in terms of the infinitesimal Lie algebra action, 
as in Equation Q. 

(4) Choose a Lie group integrator, plug in all the building blocks, and solve the 
problem. 

Examples of Lie group integrators. We list some important Lie group inte- 
grators: 

Lie Euler: 

y„+i = exp(hf(t n , y n )) ■ y n 

Lie midpoint: 

K = hf(t n + h/2,exp(K/2)-y n ) 
y n+1 = exp(K) ■ y n 

Lie RK4: There are several similar ways of turning the classical RK4 method 
into a fourth order Lie group integrator |661 168j . The following version requires 
only two commutators: 

Ki = hf(t n ,y n ) 
K 2 = hf{t n /2,exp(K 1 /2)-y n ) 
K 3 = hf (t n + h/2, cxp(K 2 /2 - [K u K 2 }/8) ■ y n ) 
K A = hf(t n + h/2, exp(K 3 ) ■ y n ) 
y n+ i = cxp (X x /6 + K 2 /3 + K 3 /3 + K 4 /6 - [K u K 2 ]/3 - [Ki,K A ]/12) ■ y„ 
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RKMK methods: There is a general procedure to turn any classical Runge- 
Kutta method into a Lie group integrator of the same order [67j [48] . 

Let {a,ij,bj,Ci}j j =1 be the coefficients of a Runge-Kutta method of order p. 
The following method is a Lie group method of order p. 

for i = 1, s 

U i = I!] i a i,j K j 

Ki = dexp^ 1 (hf(exp(Ui)-y n )) 

end 

2/n+i = ex P h K i) -Vn » 

where 

°° R 11 
dexp^(K) = E ^ n (U)(K) =K- -[U,K] + -[U, [U,K]} + ■■■ 

n=0 

is the inverse of the right trivialized tangent map of the exponential, see |48j . 

Crouch Grossman and commutator free methods: Commutators pose a 
problem in the application of Lie group integrators to stiff problems, since the 
commutator often increases the stiffness of the equations dramatically. Crouch- 
Grossman [S7J E2] > and more generally commutator-free methods [TH] , avoid them 
by doing basic time-stepping using a composition of exponentials. An example of 
such a method is CF4 [19] : 

K\ = hf(t n ,y n ) 
K 2 = hf(t n /2,exp(K 1 /2)-y n ) 
K-i = hf (t„ + h/2, cxp{K 2 /2) ■ y n ) 
K 4 = hf(t n + ty2,exp(.Ki/2) • exp(K 3 - K 1 /2) ■ y n ) 
y n+1 = exp (K t /4 + K 2 /6 + K 3 /6 - K 4 /12) ■ 

exp {K 2 /6 + K 3 /6 + K 4 /4 - K 1 /12) ■ y n 

Magnus methods: In the case where f(t,y) = f(t) is a function of time alone, 
then Q is called an equation of Lie type. Specialized numerical methods have 
been developed for such problems [JUIIIZIE]- Explicit Magnus methods can achieve 
order 2p using only p function evaluations, and they are also easily designed to be 
time symmetric. 

Examples of group actions. Some group actions of interest when applying 
Lie group integrators: 

Rotational problems: Consider a differential equation y(t) = v(t) x y(t), where 
y, v e K 2 and ||y(0)|| = 1. Since \\y(t)\\ = 1 for all t, we can take M to be the 
surface of the unit sphere. Let G — SO (3) be the special orthogonal group, con- 
sisting of all orthogonal matrices with determinant 1. Let j(t) € G be a curve such 
that 7(0) = e. By differentiating ~/{t) T ~/{t) = e, we find that 7(0) T + 7(0) = 0, 
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thus g — so (3), the set of all skew-symmetric 3x3 matrices. The infinitesimal 
Lie algebra action is left multiplication with a skew matrix, the commutator is the 
matrix commutator, and the exponential map is the matrix exponential. Written 
in terms of the infinitesimal Lie algebra action, the differential equation becomes 
y(t) = v(y)y, and we may apply any Lie group integrator. Note that for low 
dimensional rotational problems, all basic operations can be computed fast using 
Rodrigues type formulas [IB] . 

Isospectral action: Isospectral differential equations are matrix valued equations 
where the eigenvalues are first integrals (invariants of motion). Consider M — 
R nx ™ and the action of G = SO(n) on M by similarity transforms, i.e. for a £ G 
and y £ M we define a ■ y — aya T . By differentiation of the action we find the 
infinitesimal action for V £ g = so(n) as V ■ y = Vy — yV, thus for this action (m) 
becomes 



y(t) = f(t,y)-y = f(t,y)y-yf(t,y), 
where /: R x M -> g. See [TUBS] f° r more details. 

Affine action: Let G — Gl(n) x R n be the affine linear group, consisting of all 
pairs a, b where a £ R™ xn is an invertible matrix and b £ R" is a vector. The 
affine action of G on M — R" is (a, b) ■ y = ay + b. The Lie algebra of G is 
g = flt(n) x R", i.e. g consists of all pairs (V,b) where V £ W lXn and b £ R™. 
The infinitesimal action is given as (V, b) ■ y = Vy + b. This action is useful for 
differential equations of the form y(t) — L(t)y + N(y), where L(t) is a stiff linear 
part and N is a non-stiff non-linear part. Such equations are cast in the form ^ 
by choosing f(t,y) — (L(t),N(y)). Applications of Lie group integrators to such 
problems is closely related to exponential integrators. In this case it is important 
to use a commutator-free Lie group method. 

Coadjoint action: Many problems of computational mechanics are naturally for- 
mulated as Lie-Poisson systems, evolving on coadjoint orbits of the dual of a Lie 
algebra [62] , Lie group integrators based on the coadjoint action of a Lie group 
on the dual of its Lie algebra are discussed in [33] . 

Classical integrators as Lie group integrators: The simplest of all group 
actions in our setting arises when G = M = R n . We can then use vector addition 
as group operation and group action. From the definitions we find that in this case 
g = R n , the commutator is 0, and the exponential map is the identity map from R™ 
to itself. The infinitesimal Lie algebra action becomes V ■ y = V, thus ^ reduces 
to y(t) — f(t,y), where f(t,y) £ R". We see that classical integration methods 
are special cases of Lie group integrators, and all the examples of methods above 
reduce to well-known Runge-Kutta methods. 
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3.3 Geometry meets algebra 

We will discuss how the B-series of Section [2] can be generalized to Lie-Butcher 
series for analyzing exact and approximate flows on manifolds. The basic building 
blocks of the numerical methods discussed above are commutators in the Lie alge- 
bra of frame vector fields g, flows of frozen vector fields, evaluation of frozen vector 
fields, and parallel transport of tangent vectors. The Lie algebra g defines an ab- 
solute parallelism on the Lie group, which yields a parallel transport of tangent 
vectors |2]. For a vector field F : A4 — >• TM. represented by a function / : M. g, 
parallel transport between T gy ftA and T y A4 is defined as r g f(y) = f{gy). This 
transport is independent of the path between gy and y, and hence it is a parallel 
transport induced by a flat connection. Furthermore, we will see that this con- 
nection has constant torsion. This is the geometric setting of Lie-Butcher series: 
a manifold with a connection which is flat with constant torsion, giving rise to 
a post-Lie algebra [521 ISS] of vector fields on M.. The enveloping algebra of a 
post-Lie algebra is called a D-algebra [59"1 1551155] . 

Let (G, M) be a homogeneous manifold and let g denote the Lie algebra of G, 
which can be identified with right invariant vector fields (i.e. invariant derivations 
on C°°(Ai)) via We let U(g) denote the universal enveloping algebra of q, 
which similarly can be identified with higher order right invariant derivations on 
C°°(M). Let I denote the identity operator on C°°(M). We define 

U(q) m :=C°°(M)^U(q) 

& M :=C°°(M)® mQ cU( Q ) M . 

Let {di} i be a basis for g. Then / £ g M can be written as / = J2i P ® ®i- This 
represents a function /: M — > g as f(x) = J2i P{ x )di & 0, and / also acts as a 
derivation f[g] := ^ Ji f l d i g for g e C°°(M). Similarly, higher order derivations 
on C°°{M) can be represented by 

h = Y J h I ®d I= hijk "' ® d Adk ■ ■ ■ € U{g) M , (9) 

where / = k, . . .) is a multi-index. 

The space U(g) M is equipped with two operations: frozen composition g,h>-> 
gh and covariant derivation g,h*-> g[h] defined for g,h € U(g) M , g — J2i 9 1 ®di> 
h = J^j h J ® dj as 

gh = ^g I h J ®didj (10) 
i,j 

g[h]=Y,9 I dih J ®dj. (11) 
i, J 

Note that these operations are independent of the choice of basis for g. The 
covariant derivation g[h] can be understood as a (higher order) derivative of h 
as it moves under parallel transport defined by the absolute parallelism. We will 
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frequently apply the alternative notation 

9>h:=g[h] (12) 

to emphasize the similarity between this operation and the product > in a Pre-Lie 
algebra. 

Let V(U(q) m ) C U(g) M denote the (first order) derivations, defined as 
V{U{q) m ) = {f£ U{q) m I f[gh] = (f[g])h + g(f[h}] for all g,h£ U{q) m }. 
It can be shown that 

V{U{q) m )=q m . (13) 

U(g) M with the operations of frozen composition and covariant derivation satisfies 
the following fundamental relationships: For any derivation / £ T>{U{q) m ) and 
any g,h £ U(q) m we have 

g[f]eV(U(&) M ) 

f[g[h]] = (fg)[h] + (f[g])[h}. 

Such an algebraic structure is called a D-algebra in [69] . We shall see that the 
free D-algebra is the algebra of forests of ordered trees, where /, g i-> fg is the 
concatenation of forests and f,g H> f[g] is left grafting. 



3.4 Ordered trees and D-algebras. 

The set 

or avl^.M.Y....!, 

of ordered rooted trees consists of all planar rooted trees. In other words, an or- 
dered rooted tree is a tree r together with a chosen order of the branches connected 
to each vertex of r. Unlike the set T C OT of rooted trees, we do not identify 
trees who differ in the order of their branches. 

The set of ordered words of elements from OT is denoted by OF, and also 
includes the empty word. OF is called the set of ordered forests, and we write 
OFc for the set of forests colored by C. Let N = M(OT) be the non-commutative 
polynomials over OT. The linear dual N* :— Hom(N,IR) is identified with the 
infinite combinations of forests, and we write (•, •) for the pairing making forests 
orthogonal. That is, {loi,uj2} — Suj 1iUJ2 , for all uii,u)2 £ OF. 

It is sometimes convenient to allow the trees to be decorated by a set C, often 
called the set of colors. This is done via a map from the vertices of the tree to the 
set C. We write OTc and OFc for the set of trees and forests colored by C. 

A basic operation on N is the left grafting product rx: N (g) N — >• N of [69 . It is 
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defined recursively by 

I rv io = lu 
u> rx I = 

lj rv» = B + (cj), (15) 

r rv = (t rv wi)oj2 + u>i(t rv CJ2) 

(™) rv cji = r rv (w rv Wi) — (r rv lj) rv wi, 

where r is a tree and wi, 002 are forests. 

The grafting product can be used to define the Grossman-Larson product (GL 
product) o : N(g)N — > N: 

B + (w 1 oaj 2 )=wirvB + (w 2 ), (16) 

extended by linearity. 

Concatenation and left grafting gives N the structure of a D-algebra, as defined 
in [69] (see also [58] EH [56]), where the composition > is left grafting. 

Definition 3.6. Let A be a unital associative algebra with product /, g H> fg 
and unit I, equipped with a non-associative composition > : A ® A — > A such that 
I > g = g for all g € A. Write T>(A) for the set of derivations: 

V{A) = {/ G A I f> (gh) = (/ > g)h + g(f > h) for all g,he A}. 

Then A is called a D-algebra if for any derivation / s and any g S A we 

have 

(iW eD(A) 

(ii)/ > {g >h) = (fg) >h + (f> g) > /i. 
In [BS] it was shown that the D-algebra N is the free D-algebra: 



Theorem 3.7 ([69]). The vector space N = k(OTc) is the free D-algebra over 
the set C. That is, for any D-algebra A and any map v: C — > D(A) there exists a 
unique D-algebra homomorphism J- v : N A such that T v (c) — v(c) for all c G C. 

C c ► N 



3! T v 



D(A) 



A 



A D-algebra homomorphism between two D-algebras A and B is an algebra 
morphism F : A -> B such that F(V(A)) C V(B), and F(a[b]) = F(a)[F(b)]. 

By applying this theorem to the D-algebra U (g) M of differential operators on 
a homogeneous manifold (defined in Section |3.3[) , we will construct elementary 



differentials and the Lie-Butcher series (Definition 3.14 and Definition 3.151 



Post-Lie algebras. In Section [3T3] we noted that the derivations in the D-algebra 
U(g) M of differential operators could be identified with q m . In general, the deriva- 
tions in a D-algebra form what is called a post-Lie algebra, and the D-algebra can 
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be identified with the universal enveloping algebra of its post-Lie algebra of deriva- 
tions. This point of view is developed in [56j . and is also being studied further 
in an ongoing project |29j where the operad behind post-Lie and D-algebras (also 
called post-associative algebras) is explored. Post-Lie algebras were introduced 
independently by Vallette in in a different context. 



Definition 3.8. A post-Lie algebra is a Lie algebra (j4, [•,■]) equipped with a 
non-commutative, non-associative product > : A ® A — > A satisfying: 

x > [y, z] — [x > y, z] + [y, x > z] (derivation property) (17) 

[x, y]>z = a^(x, y, z) - a > (y, x, z), (18) 

where a>(x, y, z) is the associator a > (x, y, z) = x > (y t> z) — [x > y) > z. 



Notice that relation (18 1 implies that a pre-Lie algebra (Section |2.5[ ) is a post-Lie 
algebra with vanishing bracket. In [SS] it is shown that the free Lie algebra over 
rooted trees colored by a set C, equipped with a post-Lie operation derived from 
grafting of trees, is the free post-Lie algebra, and that its universal enveloping 
algebra is the free D-algebra defined above. We will not pursue this point of view 
in our present study of Lie-Butcher series. 



3.5 Lie— Butcher series 

Analogous to the B-series of Section [2j the Lie-Butcher series can be used to 
represent flows - numerical or exact - on homogeneous manifolds. To achieve this 
one combines the concept of Lie series in free Lie algebras with ideas from the 
theory of B-series. An exposition of free Lie algebras and Lie series can be found 
in the book [74] by Reutenauer. 

The free Lie algebra FLA(A) over a set A of generators is the closure of the 
generators under commutation and linear combination. In particular, we have the 
free Lie algebra FLA(OT) over the set of ordered rooted trees. A Lie series is a 
series expansion: 

S = J2S n , (19) 

n>0 

where each homogeneous component is an element of FLA(OT), i.e. the S n 's are 
Lie polynomials. 

A Lie series of particular interest to us appears when computing the pullback 
of functions along flows of vector fields on homogeneous manifolds. Let F 6 X(M) 
be a vector field with flow &t,F, and if> '■ M — > g a function. Then 

§l F i> = Fty]. (20) 

t=o 



d 
dt 
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The Taylor expansion of <£>£ F tj} around therefore takes the form of a Lie series 

00 j.n / f)n \ 

„=o n - VOT t=o J (2i) 
= V + iFfofl + ~iW]] + *j-F[-F[-FM]] + ■ • • . 

Bell polynomials. The higher order derivatives of the pullbacks can be written 
in terms of non-commutative analogs of the classical Bell polynomials of [5J. These 
polynomials have also appeared in [55J EH [55J [55] . 

Definition 3.9. Let D = R(I) be the free associative algebra over an alphabet 
X = {di}, and let d : D — >• D denote the derivation given by d(di) — di+i- The 
non-commutative Bell polynomials B n = B n (d\, . . . , d n ) G R(Z) are defined by the 
recursion 

B = l 

B n = (d l + d)B n _ 1 , n>{). 



(22) 



The first few are: 

B =1 

B 1 =d 1 

B 2 = dj + d 2 

B 3 = dl + 2did 2 + d 2 di + d 3 

B 4 = df + M\d 2 + M x d 2 di + d 2 d\ + 3d x d 3 + d 3 d x + Zd\ + d A . 

Theorem 3.10 ([65, 58]). The derivatives of the pullback of a function iji along 
the time- dependent flow <&t,F can be written as: 

d n 

— $*^ = S„(F)[# (23) 

where B n (F t ) is the image of the Bell polynomials B n under the homomorphism 
given by di 1— > i* 1 ^ -1 ) ({i — \~)th derivative) . In particular 



d n 
dF 



*t Ft 1> = B n (Fi, . . .,F n M =: B n (Fi)[ip], (24) 
where F n+1 = d n /dt n \ t=0 F. 



This result allows us to obtain a Lie series corresponding to (21) for the case 
when F is non-autonomous |65| : 

®l F <<P = Y,Bn{FM- v (25) 
n! 

71 = 

Remark 3.11. It is well known that the classical Bell polynomials [5J can be 
defined in terms of determinants. The non-commutative Bell polynomials can be 
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defined in the same way, only now in terms of a non-commutative analog of the 
determinant: the quasideterminants of Gelfand and Retakh (|37j. see also |36j). 
For example, we have 



dot 



{3-1 



-1 

~ 1 )x 2 xi -1 



{ 3 1 2 )X2 Xl 




det 



xi -1 

2X2 Xi — 1 



(S3) X 2 Xi 



= Xi + 2xix 2 + x 2 xi + x 3 
= B 3 , 

where det denotes the quasidetcrminant, computed at the circled element. Sec 
[55] for details on the quasideterminants. 

The non-commutative partial Bell polynomials B n ^ k '■= B n ^{di 1 . . . , d n - k +i) 
are defined as the part of B n consisting of words cj of length k > 0, e.g. B43 — 
id\d 2 + 2did 2 di + d 2 d\. Thus 



B„=J2 B ^- ( 26 ) 



fc=i 



A Faa di Bruno bialgebra. The non-commutative Dynkin-Faa di Bruno bial- 
gebra V is obtained by using the algebra structure of V and defining the coproduct 
Ad as 

A P (I) =I®I 

" (27) 

k=l 

This extends to all of V by the product rule Ax>(didj) — Ax>(di)Ax>(dj). For 
example 

Ax>(di) = di ® di 
A v (d 2 ) = d\ <8> d 2 + d 2 ® di 
A v (did 2 ) = dl <g> did 2 + ^1^2 <8> d?. 
Note that the coproduct is not graded by | • | 

Lemma 3.12 (|58j). The coproduct of the partial Bell polynomials is: 

a 

Av(B n , k ) = B n,e <8> B t , k . (28) 
e=i 
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Note that B n \ = d n , so (27) is a special case of (28). Summing the partial B n ^ 



over k, we find the coproduct of the full Bell polynomials: 

n 

Av(B n )=Y,B n ,k®B k . (29) 



k=l 



Using Lemma 3.12 and the fact that B n k = for k > n, one can show that V is 
a bialgebra. 

Proposition 3.13 ([58]). T> = K(X) with the non- commutative concatenation 
product and the coproduct A-p form a bialgebra T>, which is neither commutative 
nor cocommutative. 



Lie— Butcher series. The Lie-series (21 1 can also be written as the Lie-Butcher 



series for the exact flow. In general, the Lie-Butcher series Bf(a) are constructed 
to represent flows given by yo h-> yt = 4 , *(yo) : 

Mv{t))=Bf{a)[* t ](vo). (30) 

Before giving the definition of Lie-Butcher series we define the elementary differ- 
entials of a vector field F; 

Definition 3.14. Let J 7 / : N — > U(q) m be the unique D-algebra morphism given 
by Theorem |3.7| by associating • to a vector field / : M — > Q. This is called the 
elementary differentials of the vector field /. 

Note that J 7 / : N — > U{q) m is given recursively by 

(i) ^/(I)=I 

(ii) T f (B + (u))=F f (uj){f} 

(iii) T f {u x u 2 ) = J>(wi)^>(w 2 ) 

The general Lie-Butcher series are expansions of elementary differentials indexed 
over ordered rooted forests. 

Definition 3.15. A Lie-Butcher series (LB-series) is a formal series expansion 
over U(q) m : 

B f (a) = £ hMa(u>)F f (u), (31) 
weOF 

where a : N — > K. 



It turns out [58] that the Lie series ( |21[ ) can be written as 

Kf^ = E (32) 



w£OT 



where 7 are the coefficients appearing when iteratively (left) grafting • onto •. 
This is the Lie-Butcher series for the exact flow. 
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To understand how Lie-Butcher series can be used to represent numerical flows 
we conduct a closer study of the coefficients a : N — > R, and understand them as 
characters in a certain Hopf algebra. This Hopf algebra allows us to formulate the 
concept of composition of LB-series. 



3.6 Composition of Lie— Butcher series 

We would like to understand the result of composing LB-series methods in a similar 
way as we did for B-series methods in Section [273} The basic problem is to deter- 
mine whether the method <I> resulting from composing two methods $ 2 o $ 1 -both 
given by LB-series-is another LB-series, and in that case, what its coefficients are. 
Just as there is a Hopf algebra governing composition of B-series (the BCK Hopf 



algebra discussed in Section 2.3), there is a Hopf algebra Hn behind the composi- 
tion of LB-series. This Hopf algebra was studied in [69], where its properties and 
its relation to the BCK Hopf algebra was explored. An introduction can also be 
found in [55] . This Hopf algebra is the dual of a version of Grossman and Larsons 
Hopf algebras in the case of ordered trees. 

The Hopf algebra of composition. As a vector space Hn is equal to N: 
Hn = K(OT). The product is given by shuffling: 

ILULJ = LjJ = W LU I 

(33) 

(nWl) LU (r 2 W 2 ) = T\(<j)l LU T 2 L0 2 ) + T%{tiU\ LU lo 2 ) 

where t%, t<x £ OT and u>i,ui2 € OF. The coproduct is given recursively by An (I) = 
I ® I and 

A n (wt) = ujT(g>l + A n (w)lu .(/®B+)A N (S-(r)), (34) 

where r € OT, u £ OF. Here LU • : N® 4 -> N®N denotes shuffle on the left and 
concatenation on the right: (wi (g) w 2 ) LU ■ (w 3 ® u; 4 ) = (toi LU W3) £g) (u^o-^). 

Note that the shuffle product also gives rise to the shuffle Hopf algebra H s h, 
whose coproduct is given by deconcatenation 74 : 

n-i 

A c (w) = I(g)w + w(g)I-|-^r 1 ---r i (g) r i+1 • • • r„, (35) 

i=i 

where (j = t\ ■ ■ ■ r„ . 

The set of ordered forests can be generated recursively from the empty forest I 
by a magmatic operation x : N — > N on N, given by /i x (wi, ^2) — B + (oj2). For 
a forest w, write lol and ujr for the left- and right part: lu = lot, x cjr. The above 
operations can be written recursively in terms of this operation: 

Concatenation: wl = Iuj = u>, and (u>i x 012) UJ3 = u)\ x (W2W3). 
Shuffle: wLUl = Iluoj = w, and wiLUw 2 = (^iLUw 2 l) x w i_r + ( w ilLUW2) x ^ir 
Coproduct An (I) = I® I, and A N (w) = ui ® 1 + A N (w L ) LU xA N (w fl ) 
The coproduct can also be written in terms of left admissible cuts, analogous 



to the coproduct in Hbck (Theorem 2.12 1: 
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Theorem 3.16 ([69 ). The coproduct in Hn can be written as 

ceFLAC(w) 



(36) 



where uj is a forest in OT. Here FLAC(w) consists of all left admissible cuts ofoj, 
including the full cut. 



A left admissible cut differs from the admissible cuts defined in Section 2.3 an 
elementary cut c of a tree r is a selection of edges to be removed from r, chosen 
in such a way that if an edge e is removed, then all the branches on the same level 
and to the left of e must also be removed. A cut results in a collection of trees 
concatenated together to form a forest P^{t) (the pruned part), and a remaining 
tree R c el {r), containing the root. A left admissible cut c — {c%, . . . , c„} on t is a 
collection of such elementary cuts, with the property that any path from the root 
to any vertex crosses at most one cut Cj. The pruned parts from each cut together 
form the pruned part P c {t) of the left admissible cut, where the parts coming from 
different cuts are shuffled together. We also include the full cut and the empty cut, 
which results in P c {t) = r and P c {t) = I, respectively. The cutting operation is 
extended to forests oj as follows: apply the B + operation to d to get a tree, cut 
this according to the above rules, without using cuts of edges coming out of the 
root, and, finally, remove the added root from R c (uj). 

See Table [4] for some examples of the coproduct An, and see or [55] for 
further examples and properties of Hn- 





AnM 


I 


1(81 






• 


• <g) 1 + I <g> • 






• • 


• •<g)I + «®»- 


f i® •• 




I 


I(g)I + »(g)» + 


i® I 




.: 


• 1(g) I + 2»»(g 


•+»® !+•»••+] 




:. 


J»(g>i+2(g)»- 


f •<g)»«+i<g)I« 





Table 4. Examples of the coproduct An 

The main result linking Hn to LB-series is the following, which is an analog of 
the Hairer- Wanner theorem (Theorem 2.13) for B-series: 



Theorem 3.17 ([69 ). The composition of two LB-series is again a LB-series: 

B f {a)[B f {[3)]=Bf{a*p), (37) 
where * is the convolution product in Hn ■ 
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3.7 Lie—Butcher series and flows on manifolds. 

We shall see how LB-series can be used to represent numerical flows. More details 
and examples can be found in |65 l 166 ] l72l 171] l69 l 158 ] [57] . 

Flows yo H> y(t) = ^t(yo) on the manifold M. can be represented in several 
different ways. Here are three procedures, giving rise to what can be called LB- 
series of Type 1, 2 and 3: 

(1) In terms of pullback series: Find a £ G(Hn) such that 

y(y(t))=B t (a){y )m for any £ U(q) m . (38) 

This representation is used in the analysis of Crouch-Grossman methods 
by Owren and Marthinsen |72| . In the classical setting, this is called a S- 
series [70] . 

(2) In terms of an autonomous differential equation: Find j3 £ £)(Hn) such that 
y(t) solves 

y'(t) = B h (f3)(y(t)). (39) 



This is called backward error analysis (confer Section 3. 

(3) In terms of a non-autonomous equation of Lie type (time dependent frozen 
vector field): Find 7 £ g(H s h) such that y(t) solves 

y'(t) = (~S t ( 7 )(yo)) y(t). (40) 

This representation is used in [65, 66L In the classical setting this is (close 
to) the standard definition of B-series. 

The algebraic relationships between the coefficients a, (3 and 7 in the above LB- 
series are [55] : 

j3 = aoe e is eulerian idempotent in Hn- 

a = cxp°(/3) Exponential wrt. GL-product 

7 = aoY~ 1 oD Dynkin idempotent in H s h(OT). 

a — Q(l) Q-operator in H s h(OT). 

The eulerian idempotent e in a commutative, connected and graded Hopf algebra 
H is the formal series e := log (Id), where Id is the identity endomorphism and 
* the convolution product in H. The Dynkin map D is the convolution of the 
antipode S and the grading operator Y, D — S*Y, and Y^ 1 oD is an idempotent. 
See e.g. [S5J for details. The operator Q is a rescaling of the Bell polynomials: 



Q nk (dx,...,d n _ k+1 ) = ~B nk (l\dx,---,j\dj,...)= V" k(u)cj 

M=n,#(w)=k 

-A (41) 
Qn(di, • ■ ■ , d n ) — >^ Q n ,k(di, ■ ■ ■ , d 



34 



A. Lundervold and H. Z. Munthe-Kaas 



By using these relationships one can convert between the various representations 
of flows. 



Example 3.18 (The exact solution). The exact solution of a differential equation 

V'(t) = F(y(t)) 

can be written as the solution of 

y' = F r y, y(0) = y o , 

where F t — F(y(t)) € g is the pullback of F along the time dependent flow of F. Let 
dt. 1 



F t = -§jBt{j). By Proposition 4-9] the pullback is given by Bt(Qi^iExact))[F], 



so 

Yo~f Exact = Q(jBxact)[»] jExact = Y^O B + (Q Exact)) ■ 

Note that this is reminiscent of a so-called combinatorial Dyson-Schwinger equa- 
tion [35]. Solving by iteration yields 

+ 

A formula for the LB-series for the exact solution was given in J7iJ| /. We observe 
that there cannot be any commutators of trees in this expression. Therefore, in 
LB-series of numerical integrators, commutators of trees must be zero up to the 
order of the method. 



Example 3.19 (The exponential Euler method). The exponential Euler method 
148)1 can be written as follows: 

y n+1 = exp(hf(y n ))y n , 

or, by reseating the vector field f , as 

Vn+i = exp(/(y„))y„. 

This equation can be interpreted as a pullback equation of the form $(y„ + i) = 
8(exp(«))[$]y n , so 

11 

a = cxp(«) = I + • + — — + ~— H 
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(Here the Grossman-Larson product is the same as concatenation). Note that 
exp(«) = Q(*), so the Type 3 LB-series for the Euler method is simply 

lEuler = •• 

Example 3.20 (The Lie-implicit midpoint method). The Lie-implicit midpoint 
method can be presented as: 

a = /(exp(i<r)p n ) ( 42 ) 

y n+ i = exp(er)y„ (43) 
We make the following ansatz: 

v = J2 «M W = «(•)•+« (J) 2+« (V) V+a ([.,!]) [;l]+a ft] !+■■■, (44) 



i.e. that a can be written as an infinitesimal LB-series. From Equation (42), we 
get that 

Since there are no forests in this expression, we must have ct([u>, u>']) — for 
all U),L)' € OT. If we write r = B + (ti • ■ ■ Tj) then, by combining Equation (45) 
with the ansatz, we see that coefficients of the LB-series are given recursively as 

£*(•) = \, 

<x(r) = ^a(n)---a(T j ). (46) 



Therefore, 



It 11. 

^Midpoint — • + ^2(4 



3.8 Substitution and backward error analysis for 
Lie— Butcher series 

In [57] the substitution law for LB-scrics methods was developed, culminating in 
a formula that can be used to calculate the modified vector field used in backward 
error analysis. 



The substitution law. The basic idea is as for B-series (Section 2.4): We 
consider substituting a LB-series into another LB-series, e.g. t3js f (p)(ct), and the 
questions are as before: is this a LB-series, and in that case, which one? The 
result is given in terms of the substitution law, defined using the freeness of the 



D-algebra N = E(OT) (Theorem 3.7 1 
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Definition 3.21. For any map a : C — > 

D(N) Theorem 3.7 implies that there a ^ c ^ ^ 

unique D-algebra homomorphism a* : 
N — > N such that a(c) — a * c for all 
c G C. This homomorphism is called a- 
substitution. D ( N ) N 



Theorem 3.22 ([57]). The substitution law defined in Definition 3.21 corresponds 
to the substitution of LB-series in the sense that 

B Bf ( a )(P)=B f (a*P). (47) 

Calculating the substitution law. To obtain a formula for the substitution 
law we consider the dual a* of a-substitution: 

(a*P,u) = (/3,a*(w)), (48) 

called the substitution character. The dual pairing (•, •) is the one induced by 
requiring that all forests in OF are orthogonal, and we may write (a,ui) = a(ui). 
The map a* is a character for the shuffle product [57]: a£(wi IU uj 2 ) = a*(a;i) LU 
a*(w 2 )- 

The formula for the substitution law is based on the cutting of trees as in the 
coproduct An- More specifically, it is based on the dual of grafting, called pruning: 

P„(w)= Y, (^-P C M>-R»- (49) 

c€LAC{ui) 

Here the sum is over the left admissible cuts, but as opposed to the cuts in the 
formula ( [36] ) for An, the full cut is not included. 

In [57 the following inductive formula for a' was obtained: 

Theorem 3.23 ([57]). We have 

E E «:Ki))B + K(F c ( W(2) )))a(i? c K 2) )), (50) 
(u)eAc ceLAC{u {2) ) 

if lo =/= 1 and a* (I) = L Here Ac denotes the deconcatenation coproduct. 

By using the magmatic operation x on N, this can also be written as a composition 
of operators: 

a* = fio (fi x ® I) o (a* <E> a* ® a) o (J® A N ) o A c . (51) 



Here A{^ denotes the coproduct in (36) with the full cut removed, and [i denotes 
concatenation. 

Some examples of the substitution character can be found in Table [5] Many 
more examples and details can be found in [57] . 

Remark 3.24. One would like the substitution law * for LB-series to be a con- 
volution product in a Hopf or bialgebra, analogous to the substitution of B-series 
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a* M 


I 

• 


I 

"(•)• 




•• 


a(m) 2 — 




I 


a(I)« + a(«) 2 I 




.1 


a(«I)« + a(«)a(2)«« - 


- a(«) 3 «I 


I. 


<*(!•)• + q(«)q(2)«. h 


- a(.) 3 L 



Table 5. Examples of the substitution character a' 



(Theorem 2.14). One way to achieve this is by obtaining a concrete description 
of the operations in the post-Lie operad. In that case one can follow the proce- 
dure in |12j . which, roughly, is the following: The post-Lie operad has a pre-Lie 
structure (general phenomenon for augmented operads), there is an associated Lie 
algebra structure, its universal enveloping algebra is a Hopf algebra, and its dual 
is the Hopf algebra for the substitution law. This is a project currently under 
investigation [2"5] . 



Acknowledgements 

We would like to thank Martin Bordemann, Kurusch Ebrahimi-Fard, Dominique 
Manchon and Ander Murua for many valuable remarks during discussions on the 
topics of this paper. 



References 

[1] E. Abe. Hopf Algebras. Cambridge University Press, 1980. 

[2] R. Abraham, J.E. Marsden, and T.S. Ratiu. Manifolds, tensor analysis, and appli- 
cations, volume 75. Springer, second edition, 1988. 

[3] A. A. Agrachev and R.V. Gamkrelidze. Chronological algebras and nonstationary 
vector fields. Journal of Mathematical Sciences, 17(1):1650-1675, 1981. 

[4] V.I. Arnold. Mathematical methods of classical mechanics. Springer, second edition, 
1989. 

[5] E.T. Bell. Partition polynomials. Annals of Mathematics, 29(l):38-46, 1927. 



38 A. Lundervold and H. Z. Munthe-Kaas 

[6] S. Blanes, F. Casas, J. A. Oteo, and J. Ros. The Magnus expansion and some of its 
applications. Physics Reports, 470(5-6):151-238, 2009. 

[7] C. Brouder. Runge-Kutta methods and renormalization. The European Physical 
Journal C: Particles and Fields, 12(3):521-534, 2000. 

[8] C.J. Budd and M.D. Piggott. Geometric integration and its applications. In P.G. 
Ciarlet and F. Cucker, editors, Handbook of numerical analysis: Foundations of 
Computational Mathematics, volume 11, pages 35-139. Elsevier, 2003. 

[9] J.C. Butcher. Coefficients for the study of Runge-Kutta integration processes. Jour- 
nal of the Australian Mathematical Society, 3(02): 185-201, 1963. 

[10] J.C. Butcher. An algebraic theory of integration methods. Mathematics of Compu- 
tation, 26(117) :79-106, 1972. 

[11] J.C. Butcher. Numerical Methods for Ordinary Differential Equations. John Wiley 
& Sons Inc, second edition, 2008. 

[12] D. Calaquc, K. Ebrahimi-Fard, and D. Manchon. Two interacting Hopf algebras 
of trees: A Hopf-algebraic approach to composition and substitution of B-series. 
Advances in Applied Mathematics, 47(2), 2011. 

[13] M. Calvo, A. Iserles, and A. Zanna. Runge-Kutta methods on manifolds. In D.F. 
Griffiths and G.A. Watson, editors, Numerical Analysis, A.R. Mitchell 75th Birthday 
Volume, pages 57-70. World Scientific, 1996. 

[14] M.P. Calvo, A. Iserles, and A. Zanna. Numerical solution of isospectral flows. Math- 
ematics of Computation, 66(220): 1461-1486, 1997. 

[15] M.P. Calvo, A. Murua, and J.M. Sanz-Serna. Modified equations for ODEs. Con- 
temporary Mathematics, 173:63-74, 1994. 

[16] M.P. Calvo and J.M. Sanz-Serna. Canonical B-series. Numerische Mathematik, 
67(2):161-175, 1994. 

[17] P. Cartier. A primer of Hopf algebras. In P. Cartier, B. Julia, P. Moussa, and 
P. Vanhove, editors, Frontiers in number theory, physics, and geometry, volume II, 
pages 537-615. Springer, 2007. 

[18] A. Cayley. On the theory of the analytical forms called trees. Philosophical Magazine 
Series 4, 13(85), 1857. 

[19] E. Celledoni, A. Marthinsen, and B. Owren. Commutator-free Lie group methods. 
Future Generation Computer Systems, 19(3):341-352, 2003. 

[20] E. Celledoni, R.I. McLachlan, B. Owren, and G.R.W. Quispel. Energy-preserving in- 
tegrators and the structure of B-series. Foundations of Computational Mathematics, 
10:673-693, 2010. 

[21] F. Chapoton. Rooted trees and an exponential-like series. ArXiv preprint, 
math:0209104, 2002. 

[22] F. Chapoton. A rooted-trees q-series lifting a one-parameter family of Lie idempo- 
tents. Algebra & Number Theory, 3(6) :61 1-636, 2009. 

[23] F. Chapoton and M. Livernet. Pre-Lie algebras and the rooted trees operad. Inter- 
national Mathematics Research Notices, 2001(8):395-408, 2001. 

[24] P. Chartier, E. Hairer, and G. Vilmart. A substitution law for B-series vector fields. 
Technical Report 5498, INRIA, 2005. 



Algebraic structures of numerical integration 



39 



[25] P. Chartier, E. Hairer, and G. Vilmart. Numerical integrators based on modified 
differential equations. Mathematics of Computation, 76(260) :f 94f , 2007. 

[26] A. Connes and D. Kreimer. Hopf algebras, renormalization and noncommutative 
geometry. Communications in Mathematical Physics, 1 99(f ):203-242, f998. 

[27] P.E. Crouch and R. Grossman. Numerical integration of ordinary differential equa- 
tions on manifolds. Journal of Nonlinear Science, 3(l):f-33, f993. 

[28] A. Dzhumadil'daev and C. Lofwall. Trees, free right-symmetric algebras, free 
Novikov algebras and identities. Homology, Homotopy and Applications, 4(2): 165- 
190, 2002. 

[29] K. Ebrahimi-Fard, A. Lundervold, D. Manchon, H.Z. Munthe-Kaas, and J.E. Vatne. 
On the post-Lie operad. Preprint, 2011. 

[30] K. Ebrahimi-Fard and D. Manchon. A Magnus- and Fer-type formula in dendriform 
algebras. Foundations of Computational Mathematics, 9:1-22, 2009. 

[31] K. Ebrahimi-Fard and D. Manchon. Pre-Lie Butcher series. Preprint, 2011. 

[32] K. Eng0. On the construction of geometric integrators in the RKMK class. BIT 
Numerical Mathematics, 40(1):41-61, 2000. 

[33] K. Eng0 and S. Faltinsen. Numerical integration of Lie-Poisson systems while 
preserving coadjoint orbits and energy. SIAM journal on Numerical Analysis, 
39(1):128-145, 2002. 

[34] K. Eng0 and A. Marthinsen. Modeling and solution of some mechanical problems 
on Lie groups. Multibody System Dynamics, 2(l):71-88, 1998. 

[35] L. Foissy. Faa di Bruno subalgebras of the Hopf algebra of planar trees from com- 
binatorial Dyson-Schwinger equations. Advances in Mathematics, 218(1):136-162, 
2008. 

[36] I.M. Gelfand, S. Gelfand, V.S. Retakh, and R.L. Wilson. Quasideterminants. Ad- 
vances in Mathematics, 193(1):56-141, 2005. 

[37] I.M. Gelfand and V.S. Retakh. Determinants of matrices over noncommutative 
rings. Functional Analysis and Its Applications, 25(2):91-102, 1991. 

[38] M. Gerstenhaber. The cohomology structure of an associative ring. Annals of 
Mathematics, 78 (2): 267-288, 1963. 

[39] R. Grossman and R.G. Larson. Hopf-algebraic structure of families of trees. Journal 
of Algebra, 126(1):184-210, 1989. 

[40] E. Hairer. Backward analysis of numerical integrators and symplectic methods. 
Annals of Numerical Mathematics, 1(1-4):107-132, 1994. 

[41] E. Hairer. Important aspects of geometric numerical integration. Journal of Scien- 
tific Computing, 25(1):67-81, 2005. 

[42] E. Hairer, C. Lubich, and G. Wanner. Geometric Numerical Integration. Springer, 
second edition, 2006. 

[43] E. Hairer, S.P. N0rsett, and G. Wanner. Solving ordinary differential equations I: 
Nonstiff problems. Springer, 1993. 

[44] E. Hairer and G. Wanner. On the Butcher group and general multi-value methods. 
Computing, 13(1):1-15, 1974. 

[45] M.E. Hoffman. Combinatorics of rooted trees and Hopf algebras. Transactions of 
the American Mathematical Society, 355(9):3795-3812, 2003. 



40 A. Lundervold and H. Z. Munthe-Kaas 

[46] A. Iserles. Numerical methods on (and off) manifolds. In F. Cucker and M. Shub, 
editors, Foundations of Computational Mathematics, pages 180-189, 1997. 

[47] A. Iserles, A. Marthinsen, and S.P. N0rsett. On the implementation of the method 
of Magnus series for linear differential equations. BIT Numerical Mathematics, 
39(2):281-304, 1999. 

[48] A. Iserles, H.Z. Munthe-Kaas, S.P. N0rsett, and A. Zanna. Lie-group methods. Acta 
Numemca, 9:215-365, 2000. 

[49] A. Iserles and S.P. N0rsett. On the solution of linear differential equations in Lie 
groups. Philosophical Transactions of the Royal Society A: Mathematical, Physical 
and Engineering Sciences, 357(1754):983-1019, 1999. 

[50] A. Iserles, G.R.W. Quispel, and P.S.P. Tse. B-series methods cannot be volume- 
preserving. BIT Numerical Mathematics, 47(2):351-378, 2007. 

[51] A. Iserles and A. Zanna. Qualitative numerical analysis of ordinary differential 
equations. In J. Renegar, M. Shub, and S. Smale, editors, The mathematics of 
numerical analysis: 1995 AMS-SIAM Summer Seminar in Applied Mathematics, 
Utah, volume 32, page 421. American Mathematical Society, 1996. 

[52] D. Kreimer. On the Hopf algebra structure of perturbative quantum field theories. 
Advances in Theoretical and Mathematical Physics, 2(2):303-334, 1998. 

[53] B. Leimkuhler and S. Reich. Simulating Hamiltonian dynamics. Cambridge Univ 
Press, 2004. 

[54] D. Lewis and P.J. Olver. Geometric integration algorithms on homogeneous mani- 
folds. Foundations of Computational Mathematics, 2(4):363-392, 2002. 

[55] A. Lundervold. Lie-Butcher series and geometric numerical integration on mani- 
folds. PhD thesis, University of Bergen, 2011. 

[56] A. Lundervold and H.Z. Munthe-Kaas. On pre-Lie-type algebras with torsion and 
curvature. Preprint, 2010. 

[57] A. Lundervold and H.Z. Munthe-Kaas. Backward error analysis and the substitution 
law for Lie group integrators. Submitted, 2011. ArXiv preprint math:1106.1071. 

[58] A. Lundervold and H.Z. Munthe-Kaas. Hopf algebras of formal diffeomorphisms 
and numerical integration on manifolds. Contemporary Mathematics, 539:295-324, 
2011. 

[59] W. Magnus. On the exponential solution of differential equations for a linear oper- 
ator. Communications on pure and applied mathematics, 7(4):649-673, 1954. 

[60] D. Manchon. Hopf Algebras in Renormalisation. In M. Hazewinkel, editor, Handbook 
of Algebra, volume 5, pages 365-427. North Holland, 2008. 

[61] D. Manchon. A short survey on pre-Lie algebras. Available at http : //math, 
univ-bpclermont . f r/~manchon/biblio/ESI-prelie2009 . pdf , 2009. 

[62] J.E. Marsden and T.S. Ratiu. Introduction to mechanics and symmetry: a basic ex- 
position of classical mechanical systems, volume 17. Springer Verlag, second edition, 
1999. 

[63] R.I. McLachlan and G.R.W. Quispel. Six lectures on the geometric integration of 
ODEs. In R.A. DeVore, A. Iserles, and S. Endre, editors, London Math. Soc. Lecture 
Note Series, volume 284, pages 155-210. Cambridge Univ. Press, 2001. 

[64] R.I. McLachlan and G.R.W. Quispel. Geometric integrators for ODEs. Journal of 
Physics A: Mathematical and General, 39:5251, 2006. 



Algebraic structures of numerical integration 



41 



H. Munthe-Kaas. Lie-Butcher theory for Runge-Kutta methods. BIT Numerical 
Mathematics, 35(4):572-587, 1995. 

H. Munthe-Kaas. Runge-Kutta methods on Lie groups. BIT Numerical Mathemat- 
ics, 38(1):92-111, 1998. 

H. Munthe-Kaas. High order Runge-Kutta methods on manifolds. Applied Numerical 
Mathematics, 29(1):115-127, 1999. 

H. Munthe-Kaas and B. Owren. Computations in a free Lie algebra. Philosophical 
Transactions of the Royal Society of London. Series A: Mathematical, Physical and 
Engineering Sciences, 357(1754):957, 1999. 

H. Munthe-Kaas and W. Wright. On the Hopf algebraic structure of Lie group 
integrators. Foundations of Computational Mathematics, 8(2):227-257, 2008. 

A. Murua. Formal series and numerical integrators, Part I: Systems of ODEs and 
symplectic integrators. Applied numerical mathematics, 29(2):221-251, 1999. 

B. Owren. Order conditions for commutator-free Lie group methods. Journal of 
Physics A: Mathematical and General, 39, 2006. 

B. Owren and A. Marthinsen. Runge-Kutta methods adapted to manifolds and 
based on rigid frames. BIT Numerical Mathematics, 39(1):116-142, 1999. 

R.S. Palais. A global formulation of the Lie theory of transformation groups. Mem- 
oirs of the AMS, 22, 1957. 

C. Reutenauer. Free Lie algebras. Oxford University Press, 1993. 

J.M. Sanz-Serna and M.P. Calvo. Numerical Hamiltonian problems. Chapman & 
Hall/CRC, 1994. 

R. Schimming and S.Z. Rida. Noncommutative Bell polynomials. International 
Journal of Algebra and Computation, 6:635-644, 1996. 

D. Segal. Free Left-Symmetrical Algebras and an Analogue of the Poincare-Birkhoff- 
Witt Theorem. Journal of Algebra, 164(3):750-772, 1994. 

R.W. Sharpe. Differential Geometry: Cartan's generalization of Klein's Erlangen 
program. Springer, 1997. 

M. Spivak. A Comprehensive Introduction to Differential Geometry, volume 1. Pub- 
lish or Perish, third edition, 2005. 

M.E. Sweedler. Hopf algebras. W.A. Benjamin, 1969. 

P.S.P. Tse. Geometric Numerical Integration: On the Numerical Preservation of 
Multiple Geometric Properties for Ordinary Differential Equations. PhD thesis, La 
Trobe University, 2007. 

B. Vallette. Homology of generalized partition posets. Journal of Pure and Applied 
Algebra, 208(2) :699-725, 2007. 

G. Vilmart. Etude d'integrateurs geometriques pour des equations differentielles. 
PhD thesis, Universite de Geneve, 2008. 

E. B. Vinberg. The theory of convex homogeneous cones. Transactions of the Moscow 
Mathematical Society, 12:340-403, 1963. 

R.F. Warming and B.J. Hyett. The modified equation approach to the stability and 
accuracy analysis of finite-difference methods. Journal of Computational Physics, 
14(2): 159-179, 1974. 



42 A. Lundervold and H. Z. Munthe-Kaas 

[86] J.H. Wilkinson. Error analysis of floating-point computation. Numerische Mathe- 
matik, 2(f) :31 9-340, 1960. 

[87] J. Wisdom and M. Holman. Symplectic maps for the N-body problem. The Astro- 
nomical Journal, 1 02:1 528-1 538, 1991. 



